High-energy magnon dispersion in the half-filled Hubbard model: A comparison with 

La 2 Cu04 
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We use quantum Monte Carlo methods and single-mode approximation to study the magnon 
dispersion in the 2D half-filled Hubbard and phonon-coupled Heisenberg models. We find that 
in the Hubbard model with U/t < 8, high-energy magnon dispersion is similar to those observed 
in inelastic neutron scattering experiments in La2Cu04. On the other hand, our studies of a 2D 
Heisenberg model coupled to dynamic optical bond phonons, fails to reproduce the experimental 
dispersion. These results can be interpreted as evidence for intermediate U /t and charge fluctuations 
in the cuprate materials. 



PACS: 75.40.Gb, 75.40.Mg, 75.10.Jm, 75.30.Ds 

While there is still no consensus on the microscopic 
origin of superconductivity in high-temperature copper 
oxide superconductors, it is widely believed that mag- 
netic fluctuations play an important role. This has led to 
extensive studies of magnetic fluctuations in the parent 
compounds of these materials, such as La2Cu04. The 
low-energy long-wavelength properties of these layered 
compounds are well described by the non-linear sigma 
model and the renormalized classical theoryu. However, 
short wavelength or high energy spin-fluctuations may 
also be relevant to superconductivity which only sets in 
after long-range antiferromagnetic order is lost. The de- 
scription of these excitations has raised many questions. 
Are there well-defined magnons at short wavelengths? 
Is there significant spectral weight in multimagnon ex- 
citations? Do the spinons present a better description 
of short-wavelength excitations? Is there a coexistence 
of spinons and magnons? The most quantitative theo- 
retical studies have been done for the nearest-neighbor 
Heisenberg model, where one finds magnons throughout 
the Brillouin Zone, with significantly reduced speatral 
weight along the antiferromagnetic zone boundary&ll 

On the experimental side, Raman scattering has long 
provided evidence for anomalous short-wavelength spin 
fluctuations!! The line-shape of the two-magnon Raman 
spectra has not been adequately explained. Similar re- 
sults are inferred from optical absorption spectraB. Un- 
derstanding these spectra is complicated by the fact that 
light couples directly to the charge degrees of freedom 
and thus their understanding within a spin-only picture 
requires many assumptions^. In contrast, neutron scat- 
tering provides a direct probe of the spin excitationsQ'El. 
Recently, the high energy magnon dispersion has been 
measured in La2Cu04 by inelastic neutron scattering. 
One finds magnon dispersion at short wavelengths, with 
(ir/2, 7r/2) magnon energy about 13% lower than that 
at (-7T, 0)lI. This is in contrast to results for the nearest- 
neighbor Heisenberg model, where one finds the magK 
energy at (tt, 0) to be 7-10% lower than at (n/2, 7r/2)t 
The latter results have also been_aeen experimentally-ip 
the materials Cu(DCOO) 2 .4D 2 oM and S^CuaO^Uil 



Explanation for the anomalous dispersion in La2Cu04 
has evoked considerable interest. Adding a second neigh- 
bor antiferromagnetic exchange takes the spectra in the 
opposite direction, ruling out the simplest possibility. 
One mechanism that has been suggested is ring exchange 
of four electrons in a plaquetteO. Such a term can be di- 
rectly obtained from the Hubbard modeltHl and was first 
proposed for the high-T c materials to explain the Raman 
spectra more than a decade agoEl One interesting aspect 
of the ring-exchange term is that, when treated fully it 
does not change the spin-wave dispersion within linear 
spin-wave theory. Thus its effect on the dispersion is a 
purely quantum effect. Recent numerical studies have 
focused on whether such a term can explain the experi- 
mental dispersion without destabilizing Neel orderL3. 

In a recent paper, Peres and Araujdla have studied the 
2D Hubbard model using mean-field theory and obtained 
a dispersion relation similar to that observed in La2Cu04 
for an intermediate value of the on-site interaction pa- 
rameter, U/t = 6, indicating that charge fluctuations - 
ignored in the Heisenberg model - need to be taken into 
account to explain the properties of this material. Given 
the inherently approximate nature of the mean-field the- 
ory, it is important to confirm and extend this result 
using independent approaches. We have used the deter- 
minant quantum Monte Carlo (det QMC) method and 
single mode approximation (SMA) to study the magnon 
dispersion in the 2D Hubbard model at half-filling. The 
det QMC method has been extensively used to study the 
ground state properties of the 2D Hubbard model at half- 
filling where the fermion sign problem can be avoided. 
We have studied the magnon dispersion along the mag- 
netic zone boundary as a function of the on-site interac- 
tion parameter, U ft. We find that with decreasing U, the 
zone-boundary dispersion changes sign and for U/t < 8 
it becomes similar to those observed in La2Cu04. 

Coupling to phonons is another potential source of 
anomalous dispersion in the cuprates, because of the 
rather high magnon energies. We have used the stochas- 
tic series expansion QMC to study a Heisenberg model 
coupled to optical bond phonons. We find that such a 
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FIG. 1. Convergence of magnon energy to the ground state 
value as a function of the inverse temperature for a lattice of 
size 7V=8x8 at 17=6.0. Error bars are smaller than symbol 
sizes. 

model fails to reproduce the experimental spectra. While 
suggestive, this, however, does not rule out the possibil- 
ity that a more realistic treatment of the spin-phonon 
couplings can mimic the experimental results. 

The Hubbard model in two dimensions is given by the 
Hamiltonian 

H = -* E ( c L c i,° + c L c i-) 

i i 

where c\ a (ci,o) creates (annihilates) an electron with spin 
a at lattice site i. The kinetic energy term includes a sum 
over nearest neighbors (i, j) and t is the hopping integral 
between adjacent sites. U is the on-site interaction, and 
fi is the chemical potential. We shall be dealing solely 
with the half-filled band, ie., (n^j + n-^i) — 1. With 
the interaction term written in a particle-hole symmetric 
from as above, the half-filled band corresponds to set- 
ting /i = 0. Henceforth we set t — 1 and express the 
interaction parameter U in units of t. 

In the limit of large U, the Hubbard model at half- 
filling maps on to the Heisenberg model with exchange 
parameter J = At 2 /U. The Heisenberg model in 2D is 
given by the Hamiltonian 

H = jJ2Si-S> (2) 

ij 

The spin operators are related to the electron creation 
and annihilation operators by the relation 

. 1 t . 
where f are the Pauli matrices. 



The det QMC0 used to study the 2D Hubbard model 
is a finite temperature method that is based on discretiz- 
ing the imaginary time j3 — LAt and employing Trotter 
approximation to decompose the full imaginary-time evo- 
lution. This approach treats the electron-electron inter- 
actions exactly and at half-hlling, where we focus in this 
work, is able to produce results with small statistical fluc- 
tuations. Ground state expectation values are obtained 
using sufficiently large values of (3. Like most other nu- 
merical approaches, the technique is limited to finite-size 
lattices and we have been able to study primarily two 
different lattice sizes, iV=6x6 and iV=8x8. While this is 
not sufficient to do a complete finite-size scaling analysis 
to obtain thermodynamic expectation values, previous 
studies have shown that a lattice of size 7V=8x8 is large 
enough to give reasonably good estimates of the ther- 
modynamic limit, especially when measuring quantities 
away from the antiferromagnetic wave vector. 

The Heisenberg model has been studied usi«g the 
Stochastic Series Expansion (SSE) QMC methodll The 
SSE method is also a finite-temperature QMC method 
based on the importance sampling of the diagonal ma- 
trix elements of the Taylor expansion of e~@ H . Ground 
state expectation values can be obtained by using suffi- 
ciently large values of /?, and there are no approximations 
beyond statistical errocs- With the recently developed 
"operator loop update" ta, the method has proven to be 
very efficient tool for studying several different models. 

In both approaches, we have measured the equal-time 
spin-spin correlation function (SiSj) and the associated 
static spin structure factor in the momentum space given 
by 

s ^ = ^X> q,(, ~ J) <K -»h;)(»jt -™ji)> (3) 

We also evaluate the static spin susceptibility in the mo- 
mentum space given by the Kubo integral 

X(q)= dr(S(-<i,T)S(<i,0)) (4) 
Jo 

The magnon energy is calculated using the relationEl 

w(q) = 25(q)/ X (q). (5) 

We start the discussion of the numerical results by 
studying the convergence of the calculated quantities as 
a function of the inverse temperature, (3. Figure |I] shows 
a plot of the magnon energies at the edge, (0, n), and the 
center, (7r/2,7r/2), of the zone boundary as a function of 
(3 for U = 6.0. The magnon energies are seen to con- 
verge to their T = values fairly rapidly with (3. This 
convergence at relatively small values of (3 is not surpris- 
ing given that we only deal with high-energy quantities. 
Henceforth, all the results presented are for (3 = 3.0. 

Next we present the results for the magnon ener- 
gies along the entire magnetic zone boundary for U — 
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FIG. 2. The magnon dispersion along a path similar to 
Figure 3a in ref. 9 for 6x6 and 8x8 lattices at U = 6.0 and 
[3 — 3.0. Error bars are smaller than symbol sizes. 
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FIG. 3. Magnon dispersion along the magnetic zone bound- 
ary for different values of U for a iV=8x8 lattice at (5 — 3.0. 
The Heisenberg plot corresponds to U = 10.0. Error bars are 
of the order of symbol sizes. 



6.0. This value of U is consistent with estimates of 
the effective q&,site energy pin La2Cu04 obtained from 
photoemissionc3 and opticalcd spectroscopy data. Figure 
|2| shows the variation of the magnon energy for a path 
along the Brillouin zone similar to the one considered in 
Figure 3a of ref. || for two different lattice sizes. Note 
that the data near (tt, tt) and (0, 0) are strongly temper- 
ature dependent, the latter being exactly proportional 
to temperature. However, as shown in Figure [I], the 
data along the zone boundary from (n, 0) to (tt/2, it/2) 
have reached close to their zero temperature values. The 
qualitative nature of the variation of the magnon energy 
along the path is in good agreement with the experi- 
mental results and mean-field predictions. In particular, 
the magnon energy at the center of the zone boundary, 
(7t/2, 7r/2), is found to be lower than that at the edge, 
(0, 7r), as seen in the experimental data and opposite to 
the dispersion found for the Heisenberg model. The mag- 
nitude of the deviation is estimated to be « —12% of 
the iV=8x8 system, matching closely the experimental 
observation. While this quantitative agreement with ex- 
perimental data is not conclusive- the magnon energies 
were evaluated in our calculation using the SMA, and 
the finite-size effects are fairly large - it is safe to con- 
clude that the 2D Hubbard model, with a value of the 
on-site interaction energy in the range estimated from 
photoemission and optical spectroscopy data, can repro- 
duce the correct magnon dispersion along the magnetic 
zone boundary. 

We now discuss the variation of the dispersion along 
the zone boundary in the 2D Hubbard model as a func- 
tion of the on-site interaction parameter, U. Figure || 
shows the results of our simulation for a N=8x8 lattice 
at (3 = 3.0. We have focused primarily on the deviation 
of the magnon energy at the center of the zone boundary 
relative to the edges. To that end, we have scaled the 
magnon energies by their value at (0,7r). The deviation 
is seen to be maximum for the smallest value of U studied 



and decreases with increasing U - becoming essentially 
flat within statistical error at U = 8.0. Unfortunately, for 
U > 8, the results of our simulation for the 2D Hubbard 
model become too noisy. However, in this limit of large 
U, the Hubbard model maps on to the Heisenberg model. 
Hence the dispersion for larger U values is expected to be 
be qualitatively similar to that for the Heisenberg model. 
We have also shown in the figure the dispersion obtained 
for a Heisenberg model with exchange parameter J cor- 
responding to U = 10.0. The point U = 8.0 marks a 
transition to "Heisenberg-like" behavior. This value of 
the transition point is consistent with that found from 
specific heat measurements of the 2D Hubbard modelEi 
To investigate other possible sources for the nature 
of magnon energy dispersion along the magnetic zone 
boundary, we have studied the effects of spin-phonon cou- 
pling within the framework of the 2D Heisenberg model. 
The model involves coupling of the spins to dynamic op- 
tical bond phonons and is given by the Hamiltonian 

H = -J^(l + a(at +afc ))S i .Sj +LU J2 a l a b ( 6 ) 

(U) b 

where b denotes the bond between the lattice sites i and 
j and the operator a\{a,b) creates (annihilates) a phonon 
on bond b. loq represents the bare phonon frequency and 
a measures the strength of the spin-phonon interaction. 

The spin-phonon model was studied using the SSE 
QMC method. Our results indicate that over any rea- 
sonable range of spin-phonon coupling strength and bare 
phonon frequency, the magnon dispersion along the mag- 
netic zone boundary is qualitatively the same as the pure 
Heisenberg model. Figure [| shows the variation of the 
magnon energy along the same path in momentum space 
as considered in Figure | for the pure Heisenberg model 
and the spin-phonon model with a = 0.1 and ujq = 0.25 
for a lattice of size iV=8x8 at 0=5.0. The value of j3 
is chosen such that j3J is the same as (3t for the data 
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FIG. 4. Magnon dispersion along the same path as in Fig- 
ure [] for a iV=8x8 lattice at f) = 5.0 for the pure Heisen- 
berg model and the spin-phonon model with a = 0.1 and 
ujo = 0.25. 



presented in Figure | with J = tt 2 /U. This allows for 
a direct comparison between the two sets of data. Once 
again we found that while the data near (jr, n) and (0,0) 
are strongly temperature dependent, the magnon ener- 
gies along the magnetic zone boundary from (tt, 0) to 
(tt/2, 7t/2) have converged close to their ground state val- 
ues. The data suggests that the deviation in the magnon 
dispersion observed for La 2 Cu04 from that found in the 
Heisenberg model cannot be explained by the type of 
spin-phonon coupling considered here. It remains to be 
seen whether a more realistic treatment of the electron- 
phonon coupling in the cuprates can account for this be- 
havior. 

To summarize, we have used the det QMC to show that 
the experimentally observed magnon dispersion along the 
magnetic zone boundary in La2Cu04 can be reproduced 
by the 2D Hubbard model, using reasonable U/t values. 
The deviation of the magnon energy at the center of the 
zone boundary is maximum for the smallest value of U 
considered, where the charge fluctuations are the largest, 
and decreases with increasing U. For U ks 8t, the disper- 
sion is flat within statistical errors. The dispersion for 
larger values of U is qualitatively similar to that of the 
pure Heisenberg model. Thus, at U rs 8t, there is a tran- 
sition to "Heisenberg-like" behavior where the effects of 
the charge fluctuations become negligible. An interesting 
question is whether this intermediate U/t regime of the 
Hubbard model is equivalent to a Heisenberg plus ring 
exchange term. At really small U/t the electrons will be 
strongly delocalized, so that a description in terms of a 
spin Hamiltonian would break down. However, at the 
intermediate U/t, found in our calculations, a higher or- 
der perturbative treatment maybe justified leading to a 
Heisenberg plus ring exchange Hamiltonian. 
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